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Abstract. The kinetic excitation of ideal magnetohydrodynamic (MHD) 
discrete Alfven eigenmodes in the second MHD ballooning stable domain is 
studied in the presence of a thermal ion temperature gradient (ITG), using 
linear gyrokinctic particle-in-ccll simulations of a local flux tube in shifted-circle 
tokamak geometry. The instabilities are identified as a-induced toroidal Alfven 
eigenmodes (aTAE); that is, bound states trapped between pressure-gradicnt- 
induced potential barriers of the Schrodinger equation for shear Alfven waves. 
Using numerical tools, we examine in detail the effect of kinetic thermal ion 
compression on aTAEs; both non-resonant coupling to ion sound waves and wave- 
particle resonances. It is shown that the Alfvenic ITG instability thresholds (e.g., 
the critical temperature gradient) are determined by two resonant absorption 
mechanisms: Landau damping and continuum damping. The numerical results 
are interpreted on the basis of a theoretical framework previously derived from 
a variational formulation. The present analysis of properties and structures of 
Alfvenic fluctuations in the presence of steep pressure gradients applies for both 
positive or negative magnetic shear and can serve as an interpretative framework 
for experimental observations in (future) high-performance fusion plasmas of 
reactor relevance. 



1. Introduction 

Ideal magnetohydrodynamic (MHD) theory for a tokamak plasma predicts local 
stability with respect to MHD ballooning modes around flux surfaces with positive 
magnetic shear when the pressure gradient exceeds a certain threshold and in the 
entire domain of negative magnetic shear [l[ . The prospect of utilizing this property 
to achieve high plasma pressures in devices designed to create thermonuclear fusion 
conditions has raised the question how this so-called second MHD ballooning stable 
domain is affected by kinetic effects, such as wave-particle interactions. Studies 
focusing on the effect of thermal ion temperature gradient (ITG) have shown differing 
results; some workers found a destabilization of virtually the entire ideal MHD stable 
domain 0, 0, 0, 0, E| , others only a kinetic broadening of the ideal MHD ballooning 

$ Present address: Associazione EURATOM-ENEA sulla Fusione, CP 65-00044 Frascati, Rome, 
Italy. 
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unstable domain Q ■ A destabilization was also observed in the presence of an electron 
temperature gradient (ETG) @, HI- These numerical gyrokinetic studies focused 
primarily on the question of stability or instability, but comparatively little insight 
was obtained with regard to the origin of the cigcnfrcqucncies, which is important 
to understand the nature of the underlying eigenmodes. In another recent study 
it was demonstrated that, in the second MHD stable domain, discrete ideal MHD 
Alfven eigenmodes exist [9(, which may be resonantly excited, for instance, through 
interactions with a sparse energetic ion population [lOlL One may now ask the question 
whether the ITG instabilities reported in 0, H, 0,MB, 01 are related to the discrete 
Alfven eigenmodes described in [ij ■ 

These ideal MHD eigenmodes, which are known as a-induced toroidal Alfven 
eigenmodes (aTAE), are solutions of the Schrodinger equation for shear Alfven waves 
(SAW), which are trapped between potential barriers induced by the pressure gradient. 
There is an infinite number of potential barriers and, thus, an infinite number of 
branches, which are denoted as aTAE(j,p). The label j > 1 identifies the potential 
well where the dominant component is trapped, and p > is the energy quantum 
number. In those parameter regimes where the cigcnfrcqucncy of a given aTAE branch 
is comparable to the diamagnctic frequency and about three times the thermal 
ion transit frequency 3 x wti@ coupling to ion sound waves and resonant interaction 
with thermal ions is expected to have a significant effect. 

The purpose of the present paper is to examine in detail such interactions between 
aTAEs and thermal ions in the presence of an ITG. The topics dealt with in the present 
paper may be summarized as follows: 

(i) New interpretation of ITG-driven instabilities reported in [1,0] m terms of aTAEs 
@, E3]. (based on results for frequencies and mode structures obtained in the 
companion paper [ill ]) 

(ii) What determines which branch of aTAEs dominates? (highlight role of minimized 
field line bending) 

(iii) Which factors determine the mode frequency in the kinetic simulations? (analysis 
with respect to field line bending, diamagnetic drift and kinetic compression) 

(iv) Which parts of the particle distribution drive the instability, which ones damp it? 
(examine response of individual particles in pressure-curvature coupling term) 

(v) How can we relate our results to a theoretical framework based on a variational 
analysis? (relation between numerical results without scale separation and 
theoretical results based on a two-scale local mode structure) 

We utilize a reduced version of the linear gyrokinetic equations derived by Chen 



& Hasegawa [12(, which describe the dynamics of drift Alfven ballooning modes and 
treat electrons as a massless fluid. The gyrokinetic equations are solved with the linear 
Sf particle- in-cell (PIC) code AWECS [l3j. Following the approach chosen in related 
earlier studies of the second MHD stable domain (e.g., 0,0, El), we employ the shifted- 
circle model equilibrium, which describes tokamak-type configurations in terms of two 
parameters: the average magnetic shear s and the normalized pressure gradient a 
|14j |. An isotropic velocity distribution is used and the modulation of the magnetic 
field strength B is ignored, so B — Bq, there are no magnetically trapped particles 
and no toroidicity-induced Alfven frequency gap. However, magnetic curvature and 

§ The factor of 3 is due to the fact that the distribution function is weighted towards high energies 
by a factor £ 5 / 2 , where £ is the kinetic energy of a particle (details are given in section [5.21 1. 
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V-B drifts are properly accounted for. Parameters are adopted from 0, 0], so that we 
may revisit the cases considered there. The model and methods used are described in 
section [2J 

In this setting, a kinetic destabilization of the second MHD stable domain is seen 
as in 0, H, 0, IE M- Based on the results given in a companion paper 11 [, where the 
effect of finite Larmor radii (FLR) on aTAEs is examined in the fluid limit, we are 
able to identify all kinetic instabilities observed in the second MHD stable domain as 
oTAEs destabilized via interaction with thermal ions. This is shown in section [3] 

In section U and [51 we analyze the simulation results from the viewpoint of non- 
resonant and resonant dynamics of kinetic thermal ion compression. The discussion 
in section [6] is used to understand numerical results and explain the observations on 
the basis of a theoretical framework derived from a variational formulation, which was 
incrementally developed in 

HE EE EE EE EE HE HE H3 (see H2 for a review). In 
section [71 the results are summarized and conclusions are drawn. 



2. Physical model and numerical methods 



We employ the linear Sf gyrokinetic particle-in-cell code AWECS introduced in 
(HI, which solves a set of one-dimensional linear gyrokinetic equations based on 
the derivation by Chen & Hasegawa [l^l for a local flux tube described by the 
ballooning formalism in terms of the extended poloidal coordinate — oo < 8 < oo 
^~6|. The high-/3 tokamak-like flux tube geometry is approximated 



14, 23, 24, \2 



using the well-known shiftcd-circlc model equilibrium, which is parametrized in terms 
of the average magnetic shear s and normalized pressure gradient a (l4| . The 
electromagnetic field perturbations are described by gyrokinetic Maxwell equations in 
terms of the magnetic flux function Sip, the electrostatic potential 6<j), and the magnetic 
compression 5B\\ . The linear gyrokinetic Vlasov equation governs the evolution of the 
fluctuating part, 6f s , of the total distribution function, f s = fo s + df s . Thermal ions 
(s = i) and electrons (s = e) are assumed to be described by an (isotropic) Maxwellian 
distribution, /o s = UsqFms- Electrons are approximated by a massless fluid, and the 
contributions from energetic ions are neglected in the present work. The resulting 
equations describe linear drift Alfven ballooning modes subject to kinetic thermal ion 
compression. A detailed description of the model and the numerical methods is given 
in [1 31 ] . Our formulation makes use of the following coefficients and parameters: 
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f = l + h 2 , g = cos6 + hsm6, h = s(6 - Ok) - asinfl. 

Much of the notation is standard and definitions of all quantities can be found in a 
companion paper ll| where the (FLR MHD) fluid limit of the model is analyzed. In 
the present study, we focus on modes with constant radial envelope {Ok = 0). 
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2.1. Gyrokinetic equations 

Following standard procedures, the adiabatic and convective responses are separated 
through the substitution 
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where = — • (bx v±)/u> c o s . The quantity SG S captures the comprcssional part of 
the non-adiabatic component of the particle response (in short, kinetic compression). 
Electrons are approximated as a massless fluid, so that SG C = 0. The evolution of 
kinetic ion compression, SG- U is governed by the gyrokinetic equation 12j, 

. ei/oi 



where a 
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and the source terms are 



-Qi(5Su - iaSSsi); 



5Su = J (S(j) — dtp) + w di J (5^ + - — J\w5B\\ , 

k± 

SS 2i = -^-(deXi)JM. 
qK 



(3) 



(4) 



(5) 



The gyrokinetic field equations for isotropic /oi and written in terms of 8G{ and in 
Laplace-transformed form (dt — > — iu>, with to = lj t + ij) read (l3| 
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where (...) = J d 3 w = J d£ J dnB/\v\\ \ is the velocity space integral. The vorticity 
equation (O consists of the following terms: field line bending (FLB), inertia (with 
FLR), MHD particle compression (MPC), magnetic field compression (MFC), and 
kinetic particle compression (KPC). Both MPC and MFC arc static compression 
effects associated with toroidal curvature and finite (3, whereas KPC captures dynamic 
compression. Equation |7|) is the quasi-neutrality condition, equation (|8|) is the 
perpendicular Ampere's law, and both contain additional kinetic compression terms. 

For given parameter values, only the fastest growing mode is captured due to 
the initial value simulation approach. Frequencies and growth rates are normalized 
by the Alfven frequency ujao = vao/(qRq)i an d velocities by the Alfven velocity 
v A0 = B /^/jjtfmfnEi. 
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2.2. Reduced model for numerical analysis of shear Alfven eigenmodes 

Being designed for a wide range of parameters, AWECS solves, by default, the full 
equations ^ and dS])-©. Before carrying out detailed analyses of the results 
obtained, it is useful to reduce the complexity of the model by eliminating negligible 
contributions. 

Our ordering analyses and numerical studies show that the following 
approximations can be made in the gyrokinetic Maxwell- Vlasov equations ((3J and 
©-© at the expense of only quantitative changes when shear Alfven eigenmodes are 
of interest: 

• Let uj5B\\ = —QpSip, which effectively eliminates f2 p (high-/? correction to the 
curvature drift) from Wdi at the order considered, so that the magnetic drift is 
given by ui di = £7 K (v 2 + fiBo)/ui c0i . 

• Terms containing the factor vu {dg\i)J\ (ion perturbation to the toroidal current) 
can be ignored. 

• Although the parallel electric field, SE\\ = —dg (Sep — Sip), influences the frequency 
and growth rate quantitatively, its presence does not affect Alfven eigenmodes; 
so, for the purpose of simplicity, its contribution may be ignored. 

For convenience, the following dimensionless quantities are used: 

lo = , w*i = = , il K = = qkoiVtig{0), (9) 

^AO WAO £n ^ci^AO 

V „. 2 a2 V ti A a£ n 



«ao' ~" Vti < 2 2 q2[l + r ]i + T?(l+r le )y 
where ujao — VAa/(qRo) is the Alfven frequency. In the following, the hats will be 
neglected on all quantities except fcoi- For completeness, we retained the contributions 
due to SE\\, which may be readily dropped at a later stage by letting — > 0. The 
resulting reduced model for drift Alfven ballooning modes subject to kinetic compression 
may be compactly written in the form of a Schrodingcr equation, 

J< - V eS [9)8V s = 0. (10) 

Here, = vW, / = {k±/k<,) 2 = 1 + h 2 , and SV'j = d 2 (<5* s )/d6» 2 . The effective 
Schrodingcr potential is written as 

V eS = V + V m ,u + V m , T + K,FLR + V ki . (11) 

The ideal MHD potential, V, contributions from the inertia term, V mtU and V m ^ T , and 
the FLR correction for the "MPC+MFC" term, V^flr.; were introduced and discussed 
in ll|. The contribution from the KPC term in ([6]) is captured by which is the 



only addition made here to the FLR MHD equation (7) of [ll[ : 

ktfVf w ao OOi 

The evolution of SGi = SGi x m{T- x / '(eniow ao) is governed by the reduced gyrokinetic 
equation 

[v {l dg - i(w - Wdi)] SGi = - u)F m J {u^ B + SU B )/y/f, (13) 

and the contribution of the parallel electric field, represented by 6U S = uiyff(5(f> — Sip), 
is given by the quasi-ncutrality condition in the form 

(1 + t^)SU s = - r Ti) - lu(1 - T )} SV S - r?^f(J SG)-,. (14) 
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Table 1. Default physical parameters in the two cases considered in this paper: 
one with lower magnetic shear (s = 0.4) and one with higher shear (s = 1.0). 
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The linear time-independent Schrodinger equation (|10[) is used in section 3] to define 
quadratic forms, and in section [5] to make a connection to a theoretical framework 
based on a variational formulation. 

2.3. Parameter values 

In the present paper, we analyze two cases: one with lower magnetic shear (s = 0.4) 
and one with higher magnetic shear (s = 1.0). The default parameters are listed in 
table [T] They were adopted (with minor changes) from two earlier linear gyrokinetic 
studies on Alfvenic ITG-driven instabilities [2|, |7[ . In the case adopted from Hirose 
et al. Q, fcoi is increased from 0.1 to 0.2 in order to demonstrate the excitation of 
the aTAE(l,0) branch besides the (2,0) branch. The normalized wavenumber used 
by Dong et al. Q is adjusted by a factor \/2 due to a different definition of vw 
(fc i = 0.212 = 0.3/V2). 

Note that is is not meaningful to simply compare the values of parameters (s, a) of 
the shifted-circle equilibrium model used here with parameter values in experimental 
configurations with non-circular flux surfaces. The location of the stability boundaries 
of the ideally unstable domain and the distribution of the aTAE bands in the s-a plane 
are sensitive to the flux surface geometry, so it is likely that a given region in parameter 
space of the simple model may correspond to another region of the parameter space 
in a more realistic geometry A meaningful comparison would need to consider the 
shape of the Schrodinger potential at different points in parameter space, rather than 
the values of parameters like s and a. To our knowledge, this has not been done yet, 
so we have no basis for making comparisons with experiments at this point. 

3. Gyrokinetic simulation results: ITG-driven aTAEs 

In this section, results from linear gyrokinetic simulations using AWECS [l3| are 
presented. We adopt two cases considered previously by Hirose et al. 0] and Dong 
et al. Q ■ The physical parameters are shown in table [TJ The two cases are used 
to highlight similarities and differences between results obtained for lower and higher 
values of the magnetic shear (s = 0.4 and s = 1.0, respectively). In section [37X1 we 
present results for parameter scans with respect to the normalized pressure gradient, 
a, which provides an image of the various branches of shear Alfven cigenmodes 
destabilized in the s-a plane. In section 13.21 parameter scans with respect to the 
temperature gradient and perpendicular wavenumber, r}\ and fcoi, are inspected. 

In the following, the frequencies uj r and growth rates 7 are computed by averaging 
over a suitable time interval. The associated standard deviations, Aw r and A7, are 
shown as error bars, which indicate the level of phase-space discretization noise. The 
results shown are obtained for initial perturbations with mixed parity (i.e., asymmetric 
around 9 = 0). 




Figure 1. Low shear case, s = 0.4. a-dcpendence of (a) the growth rate 7 
and (b) the frequency aj r . For comparison with the gyrokinetic (GK) simulation 
results of AWECS (O ), (a) shows the growth rate of ballooning modes (BM) in 
ideal MHD ( ) and FLR MHD ( — ■ — ), and (b) the frequencies CJBM; ^ItAE 

and hj^tae of BMs and "TAEs obtained with the FLR MHD model (— • — ) 

[ill . In addition, the diamagnetic frequency w r „i <x y/a ( ) is plotted in 

(b). The thermal transit frequency, uti, is about a factor 5 smaller than ui* p i 
[cf. equation J9j]. The dotted vertical lines separate five regions, labelled (I)-(V), 
where different branches of Alfven eigenmodes dominate. 



The numerical parameters are chosen as follows. The number of phase-space 
markers is N m = 2048 x 13 (i.e., 26624 markers distributed over 13 periods, 
6 £ [-137r,137r])l There are N g = 2048 grid points in the simulation domain 
8 G [— 0max, #max] with max = 120, and the time step for the 4th-order Runge-Kutta 
solver is At = 0.005wT^. Convergence tests were carried out to ensure acceptable 
numerical accuracy [13]. The role of boundary conditions was discussed in the 
Appendix of 11 1. 



3.1. Pressure gradient scan 

The results obtained by scanning the parameter a through the first into the second 
MHD stable domain are shown in figures [T] and [2] for the two cases in table [TJ s = 0.4 
and s = 1.0, respectively. In the domain a < 0.4 of figure [21 electrostatic ITG- 
driven modes (\8<p\ 3> \Sip\) are dominant. These were described in j2?| and will 
not be considered here. The domain where Alfvenic instabilities (\Stp\ ^> \6<f> — 5ip\) 
are observed is divided into five regions, labelled (I)-(V) in figures [T] and [2j In the 
following paragraphs, the results for regions (I)-(V) are described and discussed one- 
by-one. Results obtained with the FLR MHD model in [ll[ are used to identify the 

|| In fact, except near ideal MHD marginal stability, the results in the s = 1.0 case can be accurately 
reproduced with markers loaded in as little as 3 periods, 8 £ [— 3tt, 3tt], since oTAE eigenvalues are 
determined by the bound state components only. For the s = 0.4 case, 5 periods may be sufficient. 
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various branches of eigenmodes excited by wave-particle interactions. 

Region (I): Alfvenic instabilities are found near the first MHD ballooning stability 
boundary, 0.3 < a < 0.4 in figure[JJ and 0.5 < a < 0.6 in figure^ The frequencies are 
comparable to the diamagnetic frequency. u> t ~ o;* p i, and mode structures (not shown) 



are very extended in 9 18 1. Despite the large simulation domain used in this study 
($max = 120), the results are not numerically converged with respect to # ma x- The 
convergence test presented in figure 11 in [13| indicates that both ca r and 7 may not 
be accurate (presumably lower). Thus, we are not yet able to determine the critical 
a value for the onset of Alvenic instability below the first MHD ballooning stability 
boundary and defer a detailed study of these modes. Note that results from earlier 
studies also differ 0, US [2§| ■ This regime should be studied with a model including 
electron inertia and/or finite collisionality in order to resolve the short radial scales 
which Alfven eigenmodes near ideal MHD marginal stability exhibit (and for which 
the ion polarization current vanishes). 



Region (II): In the FLR MHD limit [llj, unstable ballooning modes (BM) arc found 
in region (II): 0.45 < a < 1.3 in figure [TJ and 0.6 < a < 2.4 in figure [2] Due to 
the stabilizing effect of FLR terms |30|, this domain is smaller than the ideal MHD 
ballooning unstable domain, which extends over the range 0.4 < a < 1.5 in figure [TJ 
and 0.6 < a < 2.6 in figure [2] 

Away from FLR MHD marginal stability, the growth rates in region (II) are 
similar for both the gyrokinetic simulation (7, O ) and FLR MHD (7, — • — ). The 
frequencies in the gyrokinetic simulation (w T , O ) are larger than in FLR MHD. 
This up-shift is a manifestation of non-resonant kinetic compression effects. The 
results in figures [TJ and figure [2] are also consistent with the expectation that the 




Figure 3. Mode structures |<5>I' S (0)| of ITG-driven aTAEs for several values of 
a. (a): s = 0.4, modes in regions (IVa) and (V) of figurc[T] (b): s = 1.0, modes 
in region (IVa) of figure [2] 



effect of kinetic thermal ion compression (both resonant and non-resonant) becomes 
increasingly important towards the FLR MHD ballooning stability boundaries where 
field line bending and ballooning-interchange drive tend to balance flcl ] . 



Region (III): Near the second FLR MHD stability boundary, ballooning instabilities 
modified by kinetic compression smoothly connect to ITG-driven aTAE(l,0). This 
branch provides the dominant instability in region (III) , which extends over the range 
1-3 S ct < 1.6 in figure [TJ and 2.4 < a < 4.6 in figure El We claim that this is an 
aTAE(l,0), because of the frequency obtained with the gyrokinetic simulation (u) T , O ) 
matches that found with the FLR MHD model (w^ae' — ' — ) [HI- This similarity 
indicates that coupling to ion sound waves has only a small effect on the frequency 
of the aTAE(l,0) branch in region (III), which is confirmed in section @] below. The 
mode structure of the instability (not shown) peaks in the central potential well (j = 1) 
and closely resembles that of the corresponding FLR MHD solution for an aTAE(l,0) 
(figure 3 in [H[). In region (III), the oTAE(l,0) branch is quasi-marginally stable 
(i.e., continuum damping is negligible) [Tl| . so we conjecture that the decrease in the 
growth rate 7(a) with increasing a is due to a reduction in the resonant drive as field 
line bending and, thus, u) r (a), increases faster than ot>* p i(a) oc i> t i(a). This detunes 
the resonance. 



Region (IVa): After a sharp drop in the frequency, the ITG-driven aTAE(2,0) branch 
becomes the dominant instability in region (IVa), which extends over the range 
!-7 < a < 2.9 in figure [lj and 4.4 < a < 7.2 in figure El We claim that this is 
an oTAE(2,0), because the frequency obtained with the gyrokinetic simulation (cj r , 

O ) matches that found with the FLR MHD model (c^tae' — ' — ) E3- 

As was explained in llj, our analysis of the FLR MHD model allows to identify 
oTAE branches only above a certain numerical threshold, a num , which is why the 
corresponding curves, w^ 2 rj!^ E (a) in figures [T] and El do not reach the low-a boundary 
of region (IVa). The fact that uj r in region (IVa) scales like a 1 / 2 shows that, in this 
regime, the oTAE(2,0) is below the threshold a^'°^ for quasi-marginal stability as 
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defined in 11 1 . This is the type of modes first seen by Hirose et al. Q . 



The mode structure of the instability, which is shown in figure [3] and peaks in the 
second potential well (j — 2), also closely resembles that of the corresponding FLR 
MHD solution for an aTAE (2,0) (figure 4 in [HI). Although the mode structures in 
figure [3] are plotted only for 9 > 0, the modes in region (IV) have mixed parity around 
9 = 0. This shows that kinetic effects preserve the degeneracy (for 9k = 0) of odd and 
even solutions of ideal MHD aTAE(j,p) branches with j > 1. 



Region (IVb) in figure^ As a increases past the boundary between regions (IVa) 
and (IVb) in figure^ a « 7.2, the frequency obtained with the gyrokinetic simulation 
(cj r , O ) continues to follow the aTAE(2,0) frequency obtained with the FLR MHD 

model {(^aTAEy — ' — )■ The mode structure also continues to resemble an aTAE(2, 0), 
as can be seen in figure Elb). Although the growth rate 7 decreases with increasing 
a, the amplitude of the large-(9 tail of the mode structure also decreases compared 
to that of the bound state component (it < \9\ < 37r). This is consistent with the 
results in figure 6 of [ll| and clearly shows that the aTAE(2,0) branch becomes quasi- 
marginally stable in region (IVb); i.e., continuum damping drops to a negligible level 
when a > a^' ^ ~ 7.2. 



Region (V) in figure^ As a increases past the boundary between regions (IV) and 
(V) in figure [U the bound state component (j,p) = (3,0) becomes dominant and 
the (2,0) component subdominant, as can be seen in figure [3]Ja). Accordingly, the 
frequency obtained in the gyrokinetic simulation (w r , O ) tends to diverge from the 

shooting code result for the FLR MHD limit (w^tae; — ' — ) which is locked onto the 
aTAE(2, 0) branch. The transition to the oTAE(3,0) branch occurs below the point 
where the (2,0) branch becomes quasi-marginally stable, so there is no region (IVb) 
here. In the parameter range scanned, it was not possible to find the oTAE(3,0) 
branch with the FLR MHD shooting code, presumably because the corresponding 
numerical threshold is too large, a num > 3.2. 



Discussion regarding continuum damping. As in the example shown in figure the 
mode structures in regions (III), (IVa) and (V) all exhibit long tails at large \9\ which 
correspond to outward propagating FLR continuum waves [ll|, a part of which gets 
reflected at the unphysical boundaries of the simulation domain (cf. Appendix of |ll|). 
The presence of this propagating component implies that the frequency of the discrete 
AlfVen eigenmode overlaps with the continuous spectrum and does not reside inside 
a gap. Consequently, in regions (IVa) and (V) of figure [T] and region (IVa) of figure 
[2] where aTAEs are below the threshold ao for quasi-marginal stability, continuum 
damping can be expected to contribute to the threshold for ITG instability. 



3.2. Dependence on temperature gradient rj^ and wavenumber fc i 

The instabilities observed in the second MHD stable domain in figures [T] and [5] all 
require that the temperature gradient, 771, and the wavenumber, fco;, exceed finite 
thresholds. In figures [4] and [5j this is demonstrated for modes of the oTAE(2,0) 
branch. In both cases, the thresholds are near r/i ~ 0.5 and feoi ~ 0.05 [panels (a) and 
(c) in both figures]. Since the values of a are chosen to be in region (IVa) of figures 
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Figure 4. Parameter scans with respect to temperature gradient r)i [(a) and (b)] 
and wavenumber fcoi [( c ) an d (d)] for an aTAE(2,0) in the case with lower shear, 
s = 0.4 [region (IVa) in figure [T]. Mode structures are shown in (e) and (f). For 
comparison with the gyrokinctic simulation results (O ), the diamagnetic drift 

frequency u> tp i ( ) is shown in (b) and (d). The arrow in (d) indicates linear 

extrapolation of u r (fcoi) towards fcgi ~~ * + - 




5 10 15 20 „ , 30 5 10 15 20 . , 30 



Figure 5. Parameter scans with respect to temperature gradient r/i [(a) and (b)] 
and wavenumber fep; [(c) and (d)] for an aTAE(2,0) in the case with higher shear, 
s = 1.0 [region (IVa) in figure 121, Mode structures are shown in (e) and (f). For 
comparison with the gyrokinetic simulation results (O ), the diamagnetic drift 

frequency oj tp ; ( ) is shown in (b) and (d). The arrow in (d) indicates linear 

extrapolation of u) r (fcrji) towards fcoi ~~ * + . 



[T] and [2] (i.e., below quasi-marginal stability in the FLR MHD limit, a < q;q 2 ' ' ) ), this 
threshold can be expected to be due both continuum damping and Landau damping. 

Corresponding results for an aTAE(l.O) in the FLR MHD quasi-marginally stable 
regime [a > ag 1,0 "*; here, region (III) of figure [2], are shown in figure [H Since 
continuum damping is negligibly small for this aTAE(l,0), the threshold is due to 
Landau damping only. Nevertheless, the instability thresholds, rji ~ 1.0 and fcoi ~ 0.1 
[panels (a) and (c)], are higher than for the modes in region (IVa) shown in figures 
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Figure 6. Parameter scans with respect to r/i [(a) and (b)] and fcgi [( c ) an< i (d)] 
for an ctTAE(l,0) in the case with higher shear, s = 1.0 [region (III) in figurc[2]. 
For comparison with the gyrokinetic simulation results ( O ) , eigenvalues obtained 

in the FLR MHD limit (□) [ill , and the diamagnetic drift frequency w* p i ( ) 

are shown in (b) and (d). Note that the MHD mode is in the quasi- marginally 
stable regime, so eigenvalues are easily obtained below the threshold for ITG 
instability using the methods of [ill ]. The mode structure is shown in figure 3(b) 
of [HI]. 



d] and We conjecture that this is due to the fact that the eigenfrequency ui T of 
the aTAE(l,0) in figure [5] is significantly larger than us* p i, while oj r ui* p i for the 
aTAE (2, 0) in figures H and This implies that the oTAE(l, 0) in figure[S]is subject 
to stronger Landau damping [cf. section [5. 1| . 

In all three cases shown in figures HHH the growth rate 7 peaks near r) v ~ (2 • • ■ 3) 
and fcj ~ (0.15 • ■ • 0.2). The subsequent decrease in 7 may be understood by noting 
that, with increasing fco;, large Larmor radii average over the wave field and, with 
increasing 771 and constant a, the thermal velocity decreases. In both cases, the drive 
is reduced [cf. section HT2"] . 

Let us now turn our attention to the frequencies in the /c i-scans [panel (d) in 
figures 0H6] . It can be seen that a naive linear extrapolation of a> r (A;oi) to fco; — * + 
yields finite frequencies in all cases. In figures an d Eld), this is indicated by 
arrows; whereas in figure Etd), we may follow the data points obtained with the FLR 
MHD model (□). Obviously, a part of this frequency is due to the finite frequency of 
ideal MHD aTAEs. In addition, as /c i decreases and the diamagnetic frequency shift 
becomes smaller, one may expect to observe the effect of coupling to ion sound waves. 
The associated frequency shift is likely to be responsible for part of the finite frequency 
obtained for fco; — ► + . This conjecture is supported by the analysis presented in 
section |4] below and by the mode structures in figure (He) and (f). 

In figure E^c), the amplitude of the continuum wave (large- \6\ tail) increases as 
77; (and, thus, 7) decreases, which indicates strong continuum damping. In contrast, 
in figure [5jf), the amplitude of the continuum wave decreases as A;oi (and both 7 
and cjj.) decreases, which indicates that the influence of continuum damping decreases 
despite the reduced ITG drive. Based on this observation, we conjecture that the mode 
frequency in[5jd) enters the kinetic thermal ion (KTI) gap [l8, 21 1, as fcoi decreases. 

Note that the situation is different in the case with lower shear shown in figure 
[U The mode structures shown in figure BJf) indicate a transition from an aTAE with 
dominant (2,0) bound state component to one with dominant (3,0) component, which 
is subject to stronger continuum damping. 
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4. Numerical analysis: SAW equation in quadratic form 

In this section, we compute components of the SAW equation in quadratic form and use 
them to obtain further insight into the a-scans in figures [1] and [2] as well as the fcoi-scans 
in figures 3] and [5j For simplicity, the calculations are carried out using the reduced 
model introduced in section [272] with = 0; i.e., SEn = 0. These simplification have 
no significant qualitative effect and are noticeable here only through minor shifts of 
the boundaries between regions (III)-(V). 

Written in quadratic form and letting = 0, equation (|10[) becomes 

= SL, SL = £L W - SWf - SW ki ; (15) 

where the components of the Lagrangian are 

SL^ = ~ I de + V m , K ) \SV S \ 2 , (16) 

jy .i—a 



SWf = i / 



\S%\ 2 + V\S^ , (17) 



SWk = -j- / ^ d8 5V* s V ki Sy s , (18) 

for a computational domain — (9 max < 9 < 9 max . The normalization constant is given 
byiV = J^g ax d8 \S^ s \ 2 /f. This choice is motivated by the discussion in the Appendix 
of [ll[ and is crucial to make parameter scans meaningful. AWECS imposes fixed 
boundary conditions, ?/>(±# max ) = 0, buffered by a region which is 0.2 x 9 max wide 
and where artificial damping is applied. Thus, the energy flux through the boundary 
is zero: $ b = ^"M^s^sl-T^ = 0. 

SL U contains the w-dependent contributions; namely, the inertia term (with FLR) 
and the FLR correction of the "MPC+MFC" term. The magnitude of the ideal MHD 
potential energy, SWf, measures the strength of FLB versus ideal MHD ballooning- 
intcrchange drive (the standard definition is SWf = (N/2)SWf). The complex number 
SWki measures the effect of kinetic thermal ion compression, both non-resonant and 
resonant contributions. 

The imaginary part of (|15|) contains no information which is of interest here 
[except that 7 oc Im(<$Wki) in the FLR MHD stable domains], so it suffices to analyze 
its real part, which we write as 

E u - SWf - Re(SW ki ) = C o7 2 ; (19) 

where the frequency-dependent term is decomposed as Re(<5L w ) = E w — C07 2 , and 
has the form 

e u = c (fi 2 - fig), fi = w r -fii, nl = nl + nl (20) 

The terms fix and fi2 are defined as follows: 

Re(SL u ) = -cu r - r Ti)/6 0i - 2fi K (l - r Ai)/&oi>* (21) 
" « ' 

2CoOi 

- (2fi K ^ pi [l - r T 2K /(l + »»)]/6oi>* +K 2 - 7 2 ) <(1 - roj/fti)^; 

V „ ' V v ' 

c n| Co 
where (X)^ = AT -1 J_^ a * d6 X\S^ S \ 2 , and the functions Tq, Ai, Tj and T2 K are 



defined in equation (13) of [11 1 



aTAE: II. Kinetic excitation with ITG 



14 



s=0.4, k =0.2 



aTAE(j,0),a=1.6, v =0.18 




Figure 7. Analysis of quadratic forms SWf (+), E w (O ) and Re(5U / i i i) (* ) m 
the case with s = 0.4. In (a), the ^-dependence is shown, with regions (II)— (V) 
corresponding to those in figure [l] In (b) and (c), fcoi~dependence is shown for 
a = 1.6 in region (III) and a = 2.3 in region (IVa), respectively. The inset boxes 
with labels (j,p) indicate the main bound state components, starting with the 
dominant one. 



The results for the case with s = 0.4 are shown in figure [7] and those for s = 1.0 
in figure [5] Note that each point in figures [7] and [5] represents a snapshots taken at a 
given time t. In order to reduce the level of phase-space discretization noise, we use a 
larger number of markers here: N m = 8192 x 13. 

The discussion in section |4~T1 focuses on SWf, and section l4~2l on Re(5Wki)- is 
shown for completeness and C07 2 is omitted because it follows from the sum of the 
other three terms [and it is very small, except in the FLR MHD ballooning unstable 
domain, region (II)] . In section 14.31 we inspect the mode structures with respect to 
coupling between different bound state components. 

4--1- Field line bending versus ballooning-interchange drive 

The qualitative form of the a-dependence of the ideal MHD potential energy, SWf (a), 
is similar in both figures[7Ja) and[5Ja). For each aTAE branch, SWf (a) has a parabola- 
like shape and the intersections of these parabolas coincide with the points where the 
dominant instability switches branches. The minima in SWf correspond to the points 
where the growth rates are largest. 

The minimum of the parabola traversing region (II) is negative, SWf < 0. In 
this region, ideal MHD ballooning drive is stronger than field line bending and the 
aTAE(l, 0) branch is replaced by unstable ideal MHD ballooning modes, modified 
by FLR and kinetic effects. The points where FLB and ballooning drive balance, 
SWf = 0, coincide with the ideal MHD ballooning stability boundaries: a cr it,i ~ 0.4 
and a cr it,2 ~ 1-5 in figure [Jj[ a) , and a cr it,i ~ 0-6 and a C rit,2 ~ 2.6 in figure [H a) . This 
implies that FLR and kinetic effects do not modify the ideal MHD mode structures, 
only the eigenvalues. In the low-a part of region (III), SWf is still negative, but small 
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Figure 8. Analysis of quadratic forms SW{ (+), E u (O ) and Re(5H / j i i) (* ) m 
the case with s = 1.0. In (a), the a-dependence is shown, with regions (II)-(IV) 
corresponding to those in figure [2] In (b) and (c), fcoi~dependence is shown for 
a = 2.9 in region (III) and a = 6.7 in region (IVa), respectively. The inset boxes 
with labels (J,p) indicate the main bound state components, starting with the 
dominant one. 



enough to enable FLR effects to stabilize ideal MHD ballooning modes. 
4-2. Kinetic compression versus ideal MHD potential energy 

In region (II) of both figures [7^a) and[8ja), we find Re(<5Wki) > 0. This confirms 
that the up-shift in the frequencies seen in figures [TJb) and [2](b) . where we compared 
gyrokinetic simulation results (O ) with FLR MHD results (wbm, — ■ — ), is due to 
kinetic compression. 

In region (III), as a is increased well beyond the second MHD stability boundary, 
Re((Wki) becomes small compared to 6Wf in both figures[7fa) and[UJa). It can even be 
seen to drop below zero; especially, when fcoi is increased as shown in figures Efb) and 
[5Jb), which indicates that the resonant contribution dominates over the non-resonant 
contribution. Thus, the eigenmodes in region (III) tend to have an "incompressible" 
response, in the sense that coupling to the ion sound branch is small; especially, for 
wavenumbers fcoi well above ITG threshold. This is consistent with figures Hfb) and 
EJb), where we saw that in region (III), the frequency from the gyrokinetic simulation 
(O ) closely matches that of the underlying FLR MHD mode (^xae> — ' — )■ 

In region (IVa), we find Rc(SW ki ) > 0. The ratio Re(<5W ki )/<Wf reaches 10% in 
figures [7ja) and 50% in figures [UJa), which implies that kinetic compression causes 
a significant up-shift in the frequency. The fco;-scans in figures [T^b) and [HIb) show 
that, for sufficiently small wavenumbers fcoii Re(#Wki) tends to rise above 5Wf. This 
confirms our earlier conjecture that non-resonant kinetic compression is responsible 
for the finite frequency obtained in figures U and by naive extrapolation towards 
fcoi -► 0+. 



aTAE: II. Kinetic excitation with ITG 



16 



Region (V) in figure [7J(a) is similar to region (IVa), and region (IVb) in figure 
[SJa) is similar to region (III). 



4-3. Coupling of multiple bound states 

In the a-scans analyzed in figures [JJJa) and[8](a), switching between different branches 
of aTAEs is observed at the boundaries of regions (III)-(V). As pointed out above, 
this is correlated with minimization of SW{ (cf. section 14. ip . with possible small 
contributions from continuum damping and wave-particle interaction. In this section, 
we inspect the transitions seen in the fcorscans in figures [JJ(b) andJS^c). The mode 
structures are not shown explicitly; instead panels (b) and (c) in figures [7] and [S] have 
inset boxes, where we list the types of bound states (J, p) seen in the corresponding 
mode structures, starting with the dominant component. 

In figure [71(b), a transition occurs around fcoi ~ 0.19, where Rc((SW / ki) switches 
sign. For fcoi < 0.19, where Rc((5W / ki) > 0, the modes excited have complex 
mode structures which consist of a dominant aTAE(2, 0) component and also strong 
components with j > 2. In addition, the mode structure exhibits a clear but smaller 
peak in the central potential well; i.e., an aTAE(l,0) component. For koi > 0.19, 
where Re(<5Wid) < 0, the aTAE (1,0) component is dominant. The slight increase 
in SWf around fcoi ~ 0.26 coincides with the location where the primary aTAE(l,0) 
couples to a (2,0) component. 

In the case of figure [7jc), the non-ideal aiTAE(2, 0) is the dominant instability in 
the entire fco; range scanned, but the mode structure also contains minor components 
with j > 2. 

The mode structures are less complex in the case with higher shear, s = 1.0, shown 
in figure[Hb) and (c). Here, fewer bound state components are coupled together, and 
an aTAE(l, 0) and an aTAE(2, 0) can be clearly identified, respectively. For fcoi ^ 0.2 
in figurelUJc), the aTAE(2,0) acquires a minor (3,0) component. 



5. Numerical analysis: wave-particle resonance 

In this section, we examine in detail the wave-particle interactions responsible for the 
excitation of the Alfvcnic instabilities described in sections [3] and [4] by analyzing the 
properties of the relevant driving term for Alfvenic ITG (AITG) modes (l9l. |28|]. 

(^JoSG), « ^2 (Ai>) jW<y Jo(\j)6Gj = ^(Aff),; (22) 
j i 
where the subscript j labels a phase space marker. We use the normalization 

£ = £/%, uj_=wj_A>ti, v\\=v\\/vti, lu = uj/ujao; (23) 

where the hats are omitted in the following. 

The method of analyzing contributions from individual markers allows to obtain 
much more detailed insight than an inspection of the imaginary component of the 
quadratic forms in equation (|15p could give. In section 15.11 a method to eliminate 
Landau damping and measure the associated damping rate is described and used. In 
section [5^1 we analyze the dependence of |Aif|j on £ , v± and v\\. The characteristic 
energy-dependence expected for AITG instabilities is reproduced numerically and the 
instability drive is confirmed to be due to interactions with the high-energy tail of the 



equilibrium distribution [31l . [3 
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Table 2. Growth rate 7 and frequency lu v of ITG-drivcn Alfvcnic instabilities in 
regions (I), (III) and (IVa) in figure [2] Results for the case where all particles are 
included ( "all E" ) arc compared to the case where only contributions from particles 
with energies satisfying £ > £ C rit are included ("£high")- £ ls normalized by v?-, 
and uj t and 7 by ljao- 



region 


(I) 


region (III) 


region (IVa) 


a = 0.5, £ c 


■it = 2.4 


a = 2.9, £ crit = 2.9 


a = 6.0, £ crit = 2.7 


all £ 


£high 


all C ^Tiigh 


all C ^high 


7: 0.016 


0.024 


0.034 0.044 


0.092 0.098 


uj r : 0.282 


0.277 


0.961 0.935 


1.222 1.191 



5.1. Landau resonance 

The quantity Q; = (wj^— in equation ([3]) originates from the phase-space gradient 
operator, Q s = ujd£ + (qRo) 2 uj~ s 1 (k± x b) ■ V, acting on the Maxwellian distribution f 0l . 
For modes satisfying 7 <g; ct> r , the boundary in velocity space at which the direction of 
Landau damping is inverted is given by the condition Qi = 0, which can be expressed 
as a condition on the energy as 

£ crit = ! + I^ i^i, (7 « a , r); (24) 

where the normalized diamagnetic frequency is given by equation ([9]). Particles with 
energies £ > £ cr it are capable of driving an instability via inverse Landau damping. 

It is possible to eliminate the effect of Landau damping by dropping the 
contributions (AH)j of particles with energies £ < £ CT a if w r 3> 7 and if at least 
one of the following conditions is satisfied: 

(a) the frequency shift caused by kinetic compression [Re(<5Wki) in (p~9|) ] is much 
smaller than the mode frequency u> r ; and/or 

(b) the particles with energies £ < £ cr jt contribute only a negligible part of that 
frequency shift. 

Assuming that this is the case (checked a posteriori), we let (AH)j = for £j < £ cr jt 
in equation (f2"2"|) . The results obtained from such a screening test are summarized in 
table [2] for three values of a, sampling regions (I), (III) and (IVa) in figure [2j 

In all cases shown in table the frequencies change only by a small amount 
(~ 3%); presumably, because at least condition (a) is satisfied, since, for the present 
parameters, the diamagnetic frequency shift is much larger than that caused by kinetic 
compression. Hence, the value of £ cr it and the resonance condition are not affected by 
the screening test and the procedure is meaningful. Consequently, the difference in 
the growth rates shown in table [2] gives a measurement of the Landau damping rate, 
~7ld, which, for the three cases inspected, lies in the range 0.006 < — 7ld ^ 0.010. 
This is somewhat lower than or comparable to the continuum damping rate found for 



aTAEs in region (IV) of figures Q] and [2] (see figure 6 in jll| ) 



5. 2. Contributions from individual particles 

The contributions of individual particles to the kinetic compression term (tod JoSG} i is 
shown in figure[9j where the quantity |A7J|j is plotted as a function of (a) the particle 
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Figure 9. Contribution of individual marker particles to the kinetic compression 
term (lj^ JoSG) i = ■ (AH)j in equation J6}. Results are shown for a = 2.9 
in region (III) of figure [2] where an ctTAE(l,0) is the dominant instability. A 
snapshot of the quantity | A/f|j is plotted for each marker particle j as a function 
of (a) the particle energy £ = (v 2 ^ + f y )/2, (b) the perpendicular velocity v±, and 
(c) the parallel velocity vn , where all velocities are normalized by n t ;. The marker 
particles are separated into two populations: £ > £ cr it ~ 2.9 (O ) and £ < £ C iit 
(+), with £ cr it given by equation 1241 . In addition, (a) shows the expected energy 
scaling \£ - £ crit |£ 5 / 2 exp(-£) ( ). 



energy £. (b) the perpendicular velocity v±, and (c) the parallel velocity -uy. For this 
analysis, we have chosen an oTAE(l,0) from region (III) in figurc[2]at a = 2.9. Similar 
results are obtained for AITG instabilities in regions (I), (IV) and (V). 

The plots constitute snapshots taken at a certain time t. Note that the marker 
particles are loaded uniformly in energy £ , so the density of markers in the plots does 
not reflect the number density of the Maxwellian-distributed physical particles. 

Energy dependence. Equations © and (f2"2"| imply that the energy dependence of 
\AH\(£) should follow the scaling \£ ^Qi/oi^di^Siil oa\£ — £ C n t \£ 5 ^ 2 exp(— £ ) which, 

indeed, gives an accurate envelope ( ) for the data points in figure [9ja). The 

location of the minimum at £ ~ 2.9 corresponds to the critical energy £ cr n given by 
equation (|24|) . so the instability drive is due to particles with £ > 2.9, while Landau 
damping is exerted by particles with £ < 2.9. 

Dependence on v±_ . Figure[^b) shows that particles with lower vj_ tend to have larger 
contributions. The peaks shift towards lower v± when fcoi is increased (not shown). 
This behavior reflects the fact that maximal wave-particle interactions require an 
optimal ration of Larmor radii to perpendicular wavelength. In the case of figurc^b), 
we have koi = 0.2 and the strongest interactions occur around v± ~ 1. Given that 
a = 2.9, s = 1.0, and the mode peak is localized around 9 € [— 7T, it], we have \JJ ~ 3. 
Thus, Larmor radii satisfying A; = k±v±/ui C Qi = ^/JkoiVj_ ~ 0.6 are favored. 

Dependence on vu . Figurc^c) shows that particles with larger i>ii tend to have larger 
contributions. The strongest response occurs around i>|| ~ 3. One reason is that the 
energies required for the instability drive through | A_ff |j are relatively large [£ ~ 5; 
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cf. figure (Ha)], so that large compensate for the low values of v± [cf. figure [S^b)] 
dictated by the constraint A; < 1 . The actual shape of the envelope for the data points 
|Ai?|j(w||) in figure (He), is determined by the resonance conditions for passing ions, 
which has the form, 

= uii [fcii — p/(qRo)] — to (not normalized); (25) 

where p = 0, ±1, ... are the side-bands generated by toroidal curvature. 

For example, for the ccTAE(l,0) used in figure the ku spectrum is roughly 
Lorentzian-shaped, centered around zero and with a half-width around 1.5/(qRo) 
(cf. figure 4 in @), the thermal transit frequency is uj ti sa 0.14 [equation ©], and 
the mode frequency is ui T w 0.9 0.05 (figure^. Thus, for particles with ~ 3, 

which show maximal response on the driving side of the distribution {£ > £ cr it), 
equation (f2"5")) becomes 0.9 w 3(fcii — p x 0.14), with k\\ = qR^k^. Hence, for p = 
these particles interact primarily with the fcy = 0.3 component of the wavenumber 
spectrum. Interactions with components k\\ > 1 may be expected to be ineffective 
since they require high-order resonances, p > 5. 

6. Discussion 

The results of the present study may be understood within a theoretical framework 
which uses a variational formulation based on a variational principle and the 
assumption that the eigenfunctions have a two-scale structure. The scale separation 
allows to write the Lagrangian in the form of a so-called generalized fishbone-like 
dispersion relation (GFLDR), which follows from asymptotic matching between the 
solution for the "ideal region" < 1) and that for the "inertial layer" (\9\ ^> 1) 
(for a review, see pH ] ) . In our specific case, this theory unifies the concepts of kinetic 
ballooning modes (KBM), beta- induced Alfven eigenmodes (BAE) and Alfvcnic ITG 
(AITG) modes, each of which captures the role of a particular physical mechanism. 

aTAEs cannot be described by a GFLDR because there is no characteristic 
separation of scales in the mode structures. The eigenvalue is determined in the region 
between the turning points of the potential well in which the eigenmode is trapped, so 
all relevant physical processes act in the same narrow region of size A6* sa 2ir. However, 
the physical mechanisms themselves are the same as in those cases where a two-scale 
structure exists, so the conceptual framework which was constructed from the analysis 
of the GFLDR is useful to explain the properties of aTAEs in the presence of FLR 
effects and kinetic compression. 

Let us write the SAW equation in quadratic form (|15p as 

- i$ b + SLu - 5Wu = 6W{ GF ^ R ih = 8Wt. (26) 

As indicated in equation (|2"6")) . the term iA(u>), which describes the physics of the 
"inertial layer" in the standard GFLDR, includes the following contributions: 

(i) The boundary term <E>b accounts for continuum damping (i.e., resonant absorption 
and phase mixing) , which occurs whenever the frequency of an eigenmode overlaps 
with the continuous spectrum in the region where the mode amplitude is finite. 

(ii) The ideal-MHD limit of <5L W ^ uj 2 yields ideal MHD inertia. 

(iii) Sli^, also contains FLR corrections associated with the inertia term (i.e., 
diamagnetic drift) and pressure-curvature coupling (i.e., FLR ballooning). These 
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FLR effects ( "KBM physics" ) were the subject of investigation in the companion 
paper Due to FLR effects, the contribution of ion inertia decays with 

increasing 9 (decreasing radial scale length), so the Alfven velocity diverges until 
electron inertia and/or collisions (not included here) take over. 

Kinetic thermal ion compression effects, captured by <5Wki, are separated into two 
categories: non-resonant and resonant: 

(iv) The non-resonant contribution of kinetic compression ( "B AE physics" ) accounts 
for coupling between the shear AlfVen wave branch and the ion sound wave branch 
via geodesic curvature. 

(v) The resonant contribution of kinetic compression represents (inverse) Landau 
damping in the presence of a temperature gradient (here, ion temperature 
gradient, "AITG physics"). 

Together with the ideal MHD potential energy, SWf (which includes the effects of 
geometry, field line bending, and free expansion energy associated with the pressure 
gradient), the boundary term <I> b (i) and the inertia term J 1 (ii) describe ideal MHD 
shear Alfven waves and eigenmodes. The non-ideal effects (iii)-(v) influence an 
ideal MHD cigenmode by altering its frequency, growth rate and mode structure. 
In addition, they introduce gaps and, hence, various accumulation points in the 
continuous SAW spectrum, which serve as birth points for various types of discrete 
Alfven eigenmodes [2lj ]. 



Normal modes. Here and in the companion paper ll|, we have seen these 
mechanisms at work on two types of ideal MHD eigenmodes: ballooning modes and 
oTAEs. Away from the ideal MHD marginal stability boundaries, these branches were 
clearly identified by comparison with FLR MHD shooting code results from (see 
figures Q] and [2] above). Near the ideal MHD marginal stability boundaries, all terms 
in equation (|26[) tend to be small and comparable, so different physics and different 
types of eigenmodes are mixed. 

As we scan through the second MHD ballooning stable domain by increasing the 
normalized pressure gradient a for constant magnetic shear s, different branches of 
oTAEs are excited in turns. The switching between branches is found to be correlated 
with the minimization of ideal MHD potential energy, SW[ [sec figures [7ja) and [8ja)] . 

Resonant drive with ITG. The broadening of the ideally unstable domain due to 
wave-particle resonances is seen as predicted theoretically ill Til . 33| and previously 



shown numerically in local @, 0, [28|, [§9[ an d global [34l,l35f simulations. Measurements 
such as those in figure [5Jc) illustrate clearly that, the maximal drive is due to particles 
in the range 2 < v\\/vti < 4. This is consistent with the analytic estimate that can 
be made by maximizing the strength of wave-particle interactions as discussed in 
* 3lL 32|. Thus, in order to correctly estimate the importance of the transit resonance 



(e.g., when deriving model equations using formal ordering arguments), the reference 
velocity should be chosen several times larger than the thermal velocity u t j . 

The analysis in [l9|, HI] led to the conclusion that ITG-drive is effective only 
when an eigenmode is not far from ideal MHD marginal stability. This statement, 
which is based on an analysis of the GFLDR (two-scale assumption) in the first 
MHD ballooning stable domain, should be rendered more precisely in the light of 
our numerical simulations, which do not impose the scale-separation constraint and 
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show a dcstabilization of oTAEs in the second MHD stable domain, far away from 
marginal stability. What one can say in general is that eigenmodcs which minimize 
SWf are most easily destabilized via ITG-drive. This can be seen by comparing figures 
Q] and [5] with figures [TJa) and[5Ja). 

Consistently with many earlier works, the resonant dcstabilization via ITG 
is found to require a finite temperature gradient, rj\ > 0(0.5), and normalized 
wavenumbers of the order fcoi = fctf^ti/wcoi > 0(0.1) (cf. section HO)) . The finite 
instability thresholds are due to resonant wave absorption by both thermal ions 
(Landau damping, 7ld) and shear Alfvcn continuous spectrum (continuum damping, 
7cd)- Our numerical analyses in section [5.11 show that both absorption rates can 
be comparable: — 7ld — 7cd ~ 0(1O -2 ). In particular, this is true in regimes 
where the dominant instabilities are below the threshold for quasi-marginal stability 
of the corresponding aTAE(j,n) branch {a < ot^ ,p \ as defined in section 3 of fill]). 
The instabilities described in [2J belong to this class [e.g., regions (IVa) and (V) 
in figure Q] . In the case of aTAE above the threshold for quasi-marginal stability 
(a > c(q' p ^), continuum damping becomes negligible [— 7cd < 0(1O~ 4 )], so the 
instability thresholds in r]i and k i are mainly due to Landau damping [e.g., regions 
(III) and (IVb) in figure [j. 



Non-resonant contribution. Kinetic compression, Re((5W / ki), generally plays an 
important role in setting the cigenfrequency when the frequency is low enough so that 
the phase velocity of the wave becomes comparable to the typical ion sound speed, in 
which case part of the wave energy is spent on compressing and decompressing the ions. 
This effect was seen in the MHD ballooning unstable domain and in regions where 
the dominant aTAE(j,p) branch is below quasi-marginal stability (a < 0.$''^); i.e., 
weakly trapped and subject to significant continuum damping when located outside 
a gap. Such weakly trapped aTAEs couple particularly strongly to the ion sound 
branch when the perpendicular wavelength is large compared to the Larmor radius 
(fcoi ;$ 0.1); i.e., near the instability threshold where diamagnctic effects are small. 
Examples were shown, where Re((5Wki) rises above the ideal MHD potential energy 
SWf as fcoi decreases [see panels (b) and (c) in figures [7] and [5] . In one instance, an 
aTAE(2, 0) is found to delve into the KTI gap, where continuum damping disappears 
[see figure E(f)]. 



FLR effects. As noted above, full FLR effects are important for AITG instabilities; 
in fact, in a self-consistent theory, FLR corrections and kinetic compression always 
appear together. In section 5.3 of ll|, it was shown that the sharp accumulation points 
of the continuous spectrum known from lowest-order FLR models turn into smooth 
transitions when FLR effects are fully retained. One consequence is that the associated 
frequency shifts depend on the eigenmode structure. In ll|, this phenomenon was 



discussed for the diamagnetic gap, but the same applies to the KTI gap (see Appendix). 



7. Conclusion 



In this work, the linear destabilization of discrete shear Alfven eigenmodes in 
tokamak geometry has been studied numerically in the presence of kinetic thermal 
ion compression and thermal ion temperature gradient (ITG). The simulations where 
carried out using the 1-D linear Sf particle- in-cell code AWECS [l3| . which solves 
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gyrokinctic equations based on a model by Chen & Hasegawa [l2| in a local flux 
tube and an s-a model equilibrium Of particular interest was the second 

magnctohydrodynamic (MHD) ballooning stable domain 

In the scenario considered, there are two fundamental types of ideal MHD shear 
Alfven cigenmodes induced by magnetic field non-uniformities due to curvature and 
a steep pressure gradient. The first is the well-known unstable ballooning mode 
(BM), for which ideal MHD potential energy is negative, SW[ < 0. The second 
type encompasses the various branches of oscillatory (5W{ > 0) a-induced toroidal 
Alfven eigenmodes These so-called aTAEs are standing waves trapped between 
pressure-gradient-induced potential barriers. In a companion paper [llj], we described 
how aTAEs are modified by finite-Larmor-radius (FLR) effects in the fluid limit 
((SH^ki = 0). In the present work, we dealt with the effects of kinetic compression, 
(tfWti ^ 0). 

Our choice of parameters allowed us to revisit two cases which were examined in 
earlier numerical studies of Alfvenic ITG (AITG) instabilities; one with lower magnetic 
shear, s = 0.4 Q, and one with higher shear, s = 1.0 Q. The results reported by Hirose 
et al. Q for the case with lower shear, s = 0.4, are reproduced and we are now able to 
explain the observations in that and subsequent studies @, 0, 0, EL H, 0] hi terms of ITG- 
driven aTAE ground states Q and on the basis of the theoretical framework provided 
by the generalized fishbone-like dispersion relation (see [2l[ for a review) . Apart from 
clarifying earlier results, this new interpretation has two important implications: 

• The eigenfrequencies of the AITG instabilities in the second MHD ballooning 
stable domain are determined by the frequency of the ideal MHD aTAE, which 
is modified by FLR effects [llj and coupling to ion sound waves (section [4} . 

• The finite instability thresholds in temperature gradient (r]i) and wavenumbcr 
(fcoi) are due to resonant absorption by both thermal ions (Landau damping) 
and shear Alfven waves (continuum damping). Both absorption rates can be 
comparable, except for quasi-marginally stable aTAEs, for which continuum 
damping is negligible. 

Qualitatively similar results are obtained in the case with higher shear, s = 1.0: 
here the second MHD stable domain is also destabilized by ITG-driven aTAEs in the 
entire parameter range scanned. Thus, our results disagree with those presented by 
Dong et al. Q who observed a destabilization only near the second MHD stability 
boundary. 

Quasi-marginally stable aTAEs (a > ao, as defined in [IH) seem to be difficult to 
excite by interaction with thermal ions: the resulting AITG instabilities have relatively 
low growth rates and are often sub-dominant. This is presumably due to their high 
frequencies w r 3> w* p i due to strong field line bending. The excitation of these modes 
in the presence of an additional population of energetic ions has been demonstrated 
with an ideal-MHD-gyrokinetic hybrid model [13] ■ Work is underway to study aTAE 
excitation with both gyrokinctic thermal and energetic ions. The destabilization of 
shear Alfven waves in the presence of an electron temperature gradient (ETG) 0, @| 
may also be revisited and explored in the light of aTAEs. 

The results of this and related earlier studies indicate that the resonant excitation 
of aTAEs may potentially destabilize a large part of the second MHD ballooning 
stable domain, including the domain of negative magnetic shear. However, since 
the properties of aTAEs crucially depend on the structure of the SAW Schrodingcr 
potential, the present results can only serve as a proof-of-principle to motivate further 
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work in this direction. Practically relevant predictions regarding the robustness of 
the second MHD stable domain with respect to kinetic instabilities require further 
study employing more realistic high-/? equilibrium models (such as in [37] ]) and, 
eventually, global nonlinear simulations. In this context, let us note that the radial 
mode structures of local aTAE solutions have a half-width extending to neighboring 
rational surfaces (nq — m = ±1), which indicates that these modes experience strong 
toroidal coupling. This suggests that, provided the modes can be excited to significant 
nonlinear amplitudes, they may contribute to the radial transports in high-/? toroidal 
devices. 

The present analysis of properties and structures of Alfvenic fluctuations in the 
presence of steep pressure gradients applies for both positive or negative magnetic 
shear and can serve as an interpretative framework for experimental observations in 
(future) high-performance fusion plasmas of reactor relevance. 
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Appendix: FLR effects on the kinetic thermal ion (KTI) gap 

In the limit |fcoiS0| 3> L equation (|10p asymptotically approaches 

= m " + 7F^**- + *»r^*., (l*oi*<?l » i); (a.i) 

{k 0i s6) 2 k 0i s6 

where 

X(w) = (cu + t5W)(w - w«)/(l + tJ), Y{uj) = {uj- LJ* pi )2qv ti (A.2) 
Equation (|A.1[) describes FLR continuum waves and is the same as that obtained for 



the FLR MHD model because <5G; only contributes a negligible term of the form 
5*f> s sm(8)/9 2 . This may be readily verified by noting that, in the limit |fcoiS0| 3> 1, 
the gyrokinetic Vlasov equation ([3]) reduces to SGi f=a (ei/mi)QifoiJoS'tp and enters 
primarily via the term (uidJo^G) i in the vorticity equation 



Equation (|A.1[) was solved and discussed in detail in the Appendix of The 



discussion of the diamagnetic frequency shift in section 5.3 of ll| illustrates how the 
accumulation points of the shear Alfvcn continuum becomes blurred when FLR effects 
are important. The same argument can now be applied in a discussion on how FLR 
effects affect the accumulation point of the kinetic thermal ion (KTI) gap. 

When FLR effects are important, the amount of energy an eigenmodc transfers to 
an outgoing continuum wave depends on its eigenfrequency, u> T , as well as on its mode 
structure; i.e., localization along 9. The continuum damping disappears only gradually 
when the mode frequency drops below the nominal value of the upper accumulation 
point of the KTI gap obtained in the low-6; limit. In fact, since terms accounting for 
kinetic compression become negligible in the large-|fc iS^| limit [equation (jA.ip j. the 
lower bound of the KTI gap width is zero in the present model. Nevertheless, kinetic 
compression ensures that any regular eigenmode (finite fcii ) experiences only reduced 
continuum damping. 
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